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Long-Range Correlations in Self-Gravitating iV-Body Systems 
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Abstract. Observed self-gravitating systems reveal often fragmented non-equilibrium structures that feature char- 
acteristic long-range correlations. However, models accounting for non-linear structure growth are not always 
consistent with observations and a better understanding of self-gravitating TV-body systems appears necessary. 
Because unstable gravitating systems are sensitive to non-gravitational perturbations we study the effect of differ- 
ent dissipative factors as well as different small and large scale boundary conditions on idealized TV-body systems. 
We find, in the interval of negative specific heat, equilibrium properties differing from theoretical predictions 
made for gravo-thermal systems, substantiating the importance of microscopic physics and the lack of consis- 
tent theoretical tools to describe self-gravitating gas. Also, in the interval of negative specific heat, yet outside 
of equilibrium, unforced systems fragment and establish transient long-range correlations. The strength of these 
correlations depends on the degree of granularity, suggesting to make the resolution of mass and force coherent. 
Finally, persistent correlations appear in model systems subject to an energy flow. 
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1. Introduction 

Most astrophysical structures result of gravitational in- 
stabilities, from large scale cosmological structures down 
to planets. Yet, among the least understood topics in 
astrophysics we find galaxy formation and star forma- 
tion, which both involve fragmentation and the nonlin- 
ear growth of structures occurring during the non-linear 
phases of gravitational instability. 

Perhaps, one of the fundamental reason why fragmen- 
tation and structure formation via gravitational insta- 
bility appears so difficult is that we lack of consistent 
theoretical tools allowing to combine gravity with gas 
physics. Indeed, too often ignored is that classical thermo- 
dynamics does not hold for gravitating systems, because 
these are non-extensive in the thermodynamical sense 
(Lan dsberg |l972||l984| ; Tsallis |l999| Plastino & Plastino 
19991 ). Actually, many other natural systems do not re- 



Presently, the only available approach to follow the 
nonlinear phases of gravitational instabilities is to carry 
out numerical simulations. Among all the existing meth- 
ods, TV-body techniques are thought to be the most effec- 
tive to simulate self-gravitating systems as well the con- 
tinuous case as the granular phases. 

Yet, despite the considerable success of these methods 
in reproducing many observed features, many fundamen- 
tal problems remain. As mentioned above, the fragmen- 
tation and structure formation is not clearly understood. 
Related to this, CDM simulations conflict with observa- 



spect the requisites of thermodynamics. Such systems of- 
ten feature interesting phenomena such as growing long 
range correlations or phase transitions. Among the symp- 
toms of a fundamental deep problem in gravitating sys- 
tems is the appearance of negative specific heat (Lynden- 
Bell & Lynden-Bell |1977| ; Lynden-Bell |1998| ), which was 
seen for a long time as a paradox in statistical mechanics, 
since negative specific heat was thought to be impossible. 
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tions at galactic scales (Moore 1999, James et al. 2001 
Bolatto et al. [2002 ), and no theory of the ISM is presently 
able to predict the conditions of star formation. Most of 
the time the star formation process relies on recipes with 
little physical constraint. 

In situations were TV-body simulations have success 
(e.g. hot stellar systems) gravitational dynamics is suf- 
ficient to account for their main global properties, ad- 
ditional microscopic physics can be neglected. But when 
gravitational instability via fragmentation involves small 
scale physics, the outcome may be strongly dependent on 
the properties of the small scale physics. In other words, in 
situations where the growth on singularities triggered by 
gravity is allowed, the chaotic nature of gravitating sys- 
tems make them sensitive to the perturbations induced by 
non-gravitational physics. 

Therefore it is important to understand the proper- 
ties of TV-body systems subjected to various perturbations. 
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For these purposes, a numerical study of perturbed, self- 
gravitating iV-body systems is carried out. 

Among the relevant perturbations we expect that 
boundary conditions at small and large scales, as well as 
dissipative factors can play a key role. In order to char- 
acterize the individual effects of perturbations, in the tra- 
dition of analytical models, one is advised to deliberately 
use simplified models. 

A study of dissipative systems is important because 
such systems may develop long-range correlations. In the 
typical ISM, radiative cooling is very effective and induces 
a temporary energy- flow le ading the system far from equi- 
librium (Dyson & Williams 1997 ). From laboratory exper- 
iments it is well known that systems outside of equilib- 
rium may spontaneously deve lop sp atio-temporal struc- 
tures (Glansdor ff fc Prigog ine 1971 ; Nicolis & Prigogine 



1977; Prigogine 1980; Melo 1994) 



A permanent energy-flow is induced when energy loss 
due to dissipation is replenished, that is, when the system 
is continuously driven, e.g., by time-dependent bound- 
ary conditions. Such systems may develop persistent long- 
range correlations. Astrophysical examples for this are, the 
growth of structures in cosmological simulations, or the 
long term persistence of filamentary structures in shearing 
flows (Toomre & Kalnajs 1991; Huber & Pfenniger 2001a; 
Wisdom & Tremaine |98§ Salo |1995| ; Pfenniger |1998| ). 
Among other things, the effect of time-dependent poten- 
tial perturbations on dissipative self-gravitating spheres is 
studied in this paper. 

In the next Section we briefly review some theoretical 
results of the thermodynamics of self-gravitating isother- 
mal spheres. The model is presented in Sect. 3 and the ap- 
plied methods to carry out the structure analysis for the 
simulated systems are explained in Sect. 4. The results 
are presented in Sect. 5-7. In Sect. 5 quasi-equilibrium 
states of iV-body models are compared with analytical 
predictions. Sect. 6 is dedicated to a study of long-range 
correlations appearing in the interval of negative specific 
heat during the collapsing transition of gravitating sys- 
tems. Finally, in Sect. 7 the evolution of systems subjected 
to an energy-flow is discussed. 

2. Gravo-thermal Statistics 

By confining a gravitational system within a sphere, we 
can compare our numerical work with previous analyt- 
ical studies using the tools of statistical mechanics and 
point out some discrepancies. Although the theory, some- 
times referred to as gravo-thermal statistics, is not fully 
consistent, since the applied Gibbs-Boltzmann entropy is 
derived under the assumption of extensivity (Taruya & 
Sakagami EOOip , it yields some important and instructive 
results. 

Before we present the model and the results let us here 
briefl y review some theoretical findings (see Padmanabhan 



1990| for an e xtend ed discussion on the topic). 

Antonov ( |l962[) and Lynden-Bell & Wood (|l968|) stud- 



isothermal spheres and found anomalous behavior when 
compared with classical thermodynamics of extensive sys- 
tems. An example of such an anomalous behavior is the 
so called gravo-thermal catastrophe that occurs when an 
isothermal gas of energy E(< 0) and mass M is released 
within a spherical box of radius greater than Antonov's 
radius, r A = 0.335GM 2 /{-E). The analytical model pre- 
dicts that for such systems there is no equilibrium to go 
and nothing can stop the collapse of the central parts. 

Antonov's and Lynden-Bell & Wood's investigations 
based on point-like particles. Since then several studies 
were carried out with modified (i.e., non point- like) par- 
ticle potentials. Hertel & Thirring (1971) modified the 
gravitational interaction potential by introducing short- 
distance repulsive forces due to quantum degeneracy and 
Aronson & Hansen (1972) investigated the behavior of 
self-gravitating hard-spheres. Finally, Follana & Laliena 
( 2000D applied softened interaction potentials. They found 
all qualitatively the same result: Unlike the model of 
Lynden-Bell & Wood, there is always an equilibrium state 
for finite size particles. However, a phase transition, sep- 
arating a high energy homogeneous phase from a low en- 
ergy collapsing phase with core-halo structure occurs in 
an energy interval with negative microcanonical specific 
heat, 



cv = 



dT 



= -P 2 



88 



< 



(1) 



where e is the total energy and Q = 1/T is the inverse 
temperature. 

There are fewer theoretical works done with the grand 
canonical ensemble (e.g. de Vega, Sanchez & Combes 



1996 ) where mass would be allowed to be exchanged with 



the environment, therefore we limit the scope of this study 
to the canonical ensemble, that is energy can be exchanged 
with the environment, but not the mass and angular mo- 
mentum. 

Note that an ensemble allowing to exchange also an- 
gular momentum would be very relevant for astrophysi- 
cal situations, but few theoretical works have been made 



on th is important aspect (Laliena |1999| ; Flicgans & Gross 

booil). 



3. Model 

A spatially isolated, spherical iV-body system is studied. 
The N particles of the system are accelerated by gravity 
and forces induced by boundary conditions and perturba- 
tions, respectively. In all, a particle can be accelerated by 
five force types: 1.) confinement 2.) energy dissipation due 
to velocity dependent friction 3.) external forcing 4.) long 
range attraction (self-gravity) and 5.) short range repul- 
sion. T he five force types are discussed in detail in Sects. 



ied theoretically the thermodynamics of self-gravitating 
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3.1. Units 

A steep potential well confines th e particles to a sphere 
of radius R w 1 (see Sect. ^2). The total mass M 
and the gravitational constant G are equal one as well, 
G = M = 1. The dimensionless energy is the energy mea- 
sured in units of GM 2 /R, that is in our model 1. Thus 
energies and other quantities with energy dimension (like 
the temperature) can be considered as dimensionless. 



3.2. Confinement 

In order to keep the model simple and to enable a com- 
parison with established theoretical results of canonical 
isothermal spheres we apply a confinement that prevents 
gravitationally unbound particles from escape and keeps 
thus the particle number constant. 

The confinement is realized through a steep potential 
well 



conf 



oc R 



16 



(2) 



where R is the distance from the center of our spatially 
isolated system. The walls of the potential well are steep, 
so that the particles are only significantly accelerated by 
the confinement when they approach R = 1. However, the 
potential walls must not be too steep, to avoid unphysical 
accelerations near the wall, due to a discrete time-step. 

Contrary, to non-spherical reflecting enclosures the ap- 
plied potential well conserves angular momentum. 

3.3. Energy Dissipation 

The effect of different dissipation schemes is studied: Local 
dissipation, global dissipation and a scheme similar to dy- 
namical friction. For convenience we call the latter here- 
after "dynamical friction" . 

3.3.1. Local Dissipation 

The local dissipation scheme represents inelastic scatter- 
ing of interstellar gas constituents. That is, friction forces 
are added, that depend on the relative velocities and po- 
sitions of the neighboring particles. 

A particle is considered as a neighbor if its distance 
is r < A e, where e is the softening length and A is a free 
parameter. The friction force has the form: 



Fi = 



A 



( r 2 +e 2)3/2 I 



< Ae and V ■ r < 



: otherwise, 



(3) 



where V is the relative velocity of two neighboring par- 
ticles, £, rj > are free parameters and i — x,y,z. The 
term V ■ r/r = V cos d < ensures a convergent flow and 
that dissipation affects only the linear momentum. As a 
consequence angular momentum is locally conserved. 
For e <C r < A e the friction force is, 



(Vcos-d^ 
bi oc ei , 



(4) 



where ej is the ith component of the unity vector. Thus 
the friction force increases with the relative velocity and 
decreases with the distance, provided that (, 77 > 0. The 
condition r < A e ensures the local nature of the energy 
dissipation. 

3.3.2. Global Dissipation 

The global dissipation depends not on the relative veloc- 
ity V, but on the absolute particle velocity v. The global 
friction force is, 

Fi = -a Vi , (5) 

where a is a free parameter. The same friction force was 
already used in the shearing box experiments of Toomre 
& Kalnajs ( |1991| ) and Huber & Pfenniger ( |2001ajj2001b| ) . 



3.3.3. Dynamical Friction 

Setting the velocity-dispersion a = 1, which is the dis- 
persion of a virialized self-gravitating system with M = 
R = 1, Chandrasekhar dynamical friction formula can be 
parameterized in good approximation by, 



Fi = r 



(4 



(6) 



where T is a free p arameter (Chandrasekhar 1945 ; Binney 
& Tremaine 1994 ). In contrast to global dissipation this 
scheme has the feature that Fi does not increase asymp- 
totically with v. There is a maximum at v w 1.3. Thus high 
speed particles in collapsed regions dissipate no longer en- 
ergy and gravitational runaway is impeded. 

3.4. Forcing Scheme 

If the dissipated energy is replenished by a forcing scheme 
the system can be subjected to a continuous energy-flow. 
Here the energy injection is due to time-dependent bound- 
ary conditions, that is, due to a perturbation potential. 
This perturbation should on the one hand be non-periodic 
and quasi-stochastic, on the other hand it should pro- 
vide approximately a regular large-scale energy injection 
in time average. A simple perturbation potential, meeting 
these conditions, has the linear form, 



$pert(a;, y,z,t) =~f[B x (t)x + B y (t)y + B z {t)z] 



(7) 



where x,y and z are Cartesian coordinates, 7 is a free 
parameter determining the strength of the perturbation, 
and 

3 

B x{t) = ^ A Xti sm(ui x ,jt + 4> Xji ) 

z=l 

3 

Bytt) = ?, Ay,i sia(w yi it + <fr y> i) (8) 
»=i 
3 



»=i 
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where t is the time. The amplitudes Aj t i and the phases 
4>j^i are arbitrary fixed constants and remain unchanged 
for all simulations (j = x,y,z). The frequencies are given 
through u>j t i = 6Qj t i, where Sljj are arbitrary, but not 
rationally dependent constants (to avoid resonances) as 
well. Aj t i, (f>j t i and Clj^ lay down the form of the potential. 
Then, the amplitude and the frequency of the perturba- 
tion are controlled by two parameters, namely 7 and 5, 
respectively (Huber 2001). 

The perturbations are a linear combination of station- 
ary waves that do not inject momentum. 

If our system represents a molecular cloud then a per- 
turbation similar to those described above can be due to 
star clusters, clouds or other massive objects passing ir- 
regularly in the vicinity. Indeed such stochastic encounters 
must be quite frequent in galactic disks and we assume 
that the average time between two encounters is, 



1 



' pert 



fpe 



> 



"7"dyn 



(9) 



where / pcr t is the mean frequency of the encounters and 
Tdyn is the dynamical time. Thus, these encounters can 
provide a continuous low frequency energy injection on 
large scales. 

Such a low frequency forcing scheme induces at first 
ordered particle motions, which are transformed in the 
course of time to random thermal motion, due to gravita- 
tional particle interaction. 

3.5. Self-Gravity, Repulsion and Force Computation 

Interaction forces include self-gravity and repulsive forces 
whose strength can be adjusted by a parameter. Repulsive 
forces may be a result of sub-resolution processes such as 
star-formation or quantum degeneracy. The particle inter- 
action potential reads, 



*p(t) = - 



Gm 



\fr 2 + ( 2 



(r 2 + e 2 ) 



(10) 



where m is the particle mass and r is here the distance 
from the particle center. On larger scales (^ e) the poten- 
tial is a Plummer potential <!>pi (Plummer 1911; Binney 
& Tremaine 1994). On small scales (^ e) the deviation 
from a Plummer potential and the strength of the repul- 
sive force is determined by the parameter £. The interac- 
tion potential becomes for instance repulsive in the range 
|r| < e^-l) if f > 1/3 (see Fig. p. If f = 0, $ p = $n 
holds. 

The interaction forces are computed on the Gravitor 
Beowulf Cluster^] at the Geneva Observatory with a paral- 
lel tree code. This code is based on the Fortran Barnes & 
Hut (1986: 1989[ ) tree algorithm, and has been efficiently 
parallelized for a Beowulf cluster. It is available on request. 

The time-integration is the leap-frog algorithm with 
uniform time-step, which ensures the conservation of the 



-100 




-150 



-200 



Fig. 1. The potential as a function of r in units of Gm/R, 
where m is the particle mass. A test particle at \r\ < 
e \/(3£ — 1) feels a repulsive force if £ > 1/3. The size of 
the softening parameter e is indicated, as well. 



symplectic structure of the conservative dynamics. The 
time-step is, 



At < 0.1 — 



(11) 



http:/ /obswww.unige.ch/~pfennige/gravitor/gravitor.html 



where er v is the velocity-dispersion of the initial state. 

The accuracy of the force computation is given through 
the tolerance parameter, which is for the studies presented 
in this paper, 6 < 0.58. 

3.6. Code Testing 

In order to test the code, the evolution of the angular mo- 
mentum is checked. We find that the angular momentum 
J of a system with 10000 particles and local dissipation 
scheme, which is initially J = 3.4 x 10 -5 in units such that 
energy is dimensionless, remains small, J < 4 x 10~ 5 . This 
is illustrated by means of two examples. Two simulations 
with equal dissipation strength are carried out. One with 
short range repulsion, £ = 2/3, and the other without, 
£ = 0. After 10 crossing times the angular momentum is 
J = 3.3 x 10~ 5 and J = 3.8 x 10~ 5 , respectively. The de- 
viation from the initial value is larger for the simulation 
without short range repulsion. This is because the dissi- 
pation leads in this case to a stronger mass concentration, 
i.e, to shorter particle distances. Indeed, the potential en- 
ergy after 10 crossing times is U s» — 12 and U « — 17 for 
the simulation with and without short range repulsion, 
respectively. 

3.7. Parameters 

The different dissipation schemes, the forcing scheme and 
the interaction potential with the short-distance repulsive 
force have several free parameters. In order to do a rea- 
sonably sized parameter study the particle number has to 
be limited to a maximum of N = 160000. However, im- 
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portant is not the absolute particle number, but the effect 
of changing it. Thus, also N is varied. 

The model parameters and their ranges are indicated 
in Table |. 



4. Correlation Analysis 

In order to check if the perturbations of the gravitational 
system induce characteristic correlations in phase-space 
the indices D and 6 of the mass-size relation and the 
vclocity-dispersion-size relation are determined. 

The index of the mass-size relation can be found via, 



D(r) 



dlogM 
dlogr 



(r) 



(12) 



where M is the sum of the masses of all particles having 
a relative distance r. The mass size relation is then 



M(r) oc r D{r) 



(13) 



where r can be considered as the scale. For a homogeneous 
mass distribution and for a hierarchical fragmented, frac- 
tal structure the index is a constant. For the latter case 
the index is not necessarily an integer and lies in the inter- 
val, Dt < D < Z)g, where Dt and Dg are the topological 
dimension and the dim ension of the embedding space, re- 
spectively (Mandelbrot |l982| ). 

The index of the velocity-dispersion-size relation is de- 
termined by, 



S(r) 



d log ex 
dlogr 



(r), 



(14) 



where a is the velocity-dispersion of all particles hav- 
ing a relative distance r. Observations of the interstel- 
lar medium suggest a constant index = i5l on scales 
0(0.1) -0(100) pc with 0.3 <, 5 L f 0.5. This i s expressed 
by Larson 's law (e.g. Larson 1981 , Scalo 1985] , Falgarone 
& Perault |l987[ Myers & Goodman |l988j ) 



a oc r SL . (15) 

In order to preclude the effect of boundary conditions 
on the scaling relations, upper and lower cutoffs have to be 
taken into account. The lower cutoff is given by the soft- 
ening length. An upper cutoff arises from the final system 
size. Thus, the scope of application is for the velocity- 
dispersion-size and the mass-size relation, e < r < 2i?9o, 
and e < r < Rqq/2, respectively, where Rgo is the radius 
of the sphere centered at the origin which contains 90% of 
the mass (see Huber & Pfenniger 2001a). 



5. Results: I. Properties of Equilibrium States in 
Numerical and Analytical Models 

Here, some iV-body gravo-thermal experiments of systems 
with weak dissipation, i.e., of systems in quasi-equilibrium 
are presented. The results are compared with theoretical 
findings. 



~-l0 

Q. 



-15 



-20 - 
-0.5 
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Fig. 2. The Follana (n = 10) and the Plummer (e = 0.05) 
potential, in units of Gm 2 /R, where m is the particle mass. 



In order to cover a range of energy with the same ex- 
periment, a convenient way is to introduce a weak global 
dissipation scheme allowing to describe a range of quasi- 
equilibrium states. Here weak means that r^is S> Td yn , 
where Tdis is the dissipation time-scale. Then, the results 
can be compared with theoretical equilibrium states. 

Follana & Laliena ( 2000| ) examined theoretically the 
thermodynamics of self-gravitating systems with softened 
potentials. They soften the Newtonian potential by keep- 
ing n terms of an expansion in spherical Bessel functions 
(hereafter such a regularized potential is called a Follana 
potential). This regularization allows the calculation of the 
thermodynamical quantities of a self-gravitating system. 
The form of their potential is similar to a Plummer poten- 
tial with a corresponding softening length. Fig. ^| shows a 
softened Follana potential with n = 10, and a Plummer 
potential with e = 0.05. 

In their theoretical work Follana & Laliena found for 
a mild enough regularization (n < 30) a phase transition 
below the critical energy e c ~ —0.335 in a region with neg- 
ative specific heat. The transition separates a high energy 
homogeneous phase from a low energy collapsed phase 
with core-halo structure. 

We want to reproduce these findings by applying a 
Plummer potential (£ = 0). Furthermore, the effect of a 
small-scale repulsive force (£ > 1/3) is studied. 

For these purposes, simulations with a weak global dis- 
sipation strength, a = 0.025, are carried out. The dissipa- 
tion time is then Tdis = 80 Td yn = 56.6 rgp, where the free 
fall time is, ts = Tdyn 

Before we discuss the results, let us briefly present 
some model properties. The initial state is a relaxed, un- 
perturbed and confined iV-body sphere with total energy 

£ = 1. 

Assigning a particle a volume of 47re 3 /3, the volume 
filling factor is, 



V s = 



(16) 
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1/T"dyn, ■ ■ ■ ,27r/r dyll 
0,...,0.7 


Perturbation frequency 
Perturbation factor 


Perturbation potential 


A 

A 

n 
C 


0, . . . ,5000 
2 e, . . . ,50 e 
0, 1, 2 
0, 1, 2 


Local dissipation factor 
Dissipation volume 
Relative distance power 
Relative velocity power 


Local dissipation 


Q 


0, . . . ,50 


Global dissipation factor 


Global dissipation 


r 


0, . . . ,30 


Dynamical friction factor 


Dynamical firiction 


e 
t 


0.0046, 0.050 
0, 1/3, 2/3, 1 


Softening length 
Repulsion strength 


Interaction potential 


N 


6400, ...,160000 


Particle number 


Mass resolution 



Table 1. Model paramaters, characterizing energy injection (perturbation potential), dissipation, interaction potential 
(self-gravity and repulsion) and mass resolution. 




50 100 150 200 250 300 350 



Fig. 3. The potential U, kinetic T and total energy e as 
a function of time for a simulation with global dissipation 
and Plummcr softening (£ = 0). The temporal evolution 
of the virial equation is indicated as well (Viriel = 1/2). 
The dashed vertical line depicts the time, when the system 
becomes fully self-gravitating. The dash-dotted and the 
dotted line mark the interval of negative specific heat. 

where N is the particle number and R is the radius of the 
system. The volume filling factor is set Vg = 1, meaning 
that the force resolution is equal to the mass resolution. 
Then, the particle number is, N = 8000. 

The applied weak dissipation strength maintains the 

fstcm approximately thermalized and virialized (see Fig. 
. Here a system of N particles is called virialized when 
the moment of inertia I is not accelerated, 

f N 

- = 2T + ^F i r i «0, (17) 

where T is the kinetic energy, is the location of the i- 
th particle and Fi is the sum of all forces acting on the 
particle (gravitation, friction, confinement). 

The evolution of the inverse temperature (3 — 1/T as 
a function of the total energy e for two simulations with 
£ = and £ = 2/3, respectively, as well as the semi- 



analytical curve calculated by Follana & Laliena are shown 
in the left panel of Fig. ||. 

The interval of negative specific heat corresponds for 
the simulation with the Plummer potential with theoret- 
ical predictions. Also, in accordance with predictions, a 
phase transition takes place, separating a high energy ho- 
mogeneous phase, from a collapsed phase in a interval with 
negative specific heat. Yet, the simulated phase transitions 
occur at higher energies. To illustrate this, the evolution 
of the Lagrangian radii k, are shown in the middle and 
the right panel of Fig. ||. 

In the high energy homogeneous phase the system is in- 
sensitive to the short-distance form of the potential. Thus 
all systems enter the interval of negative specific heat at 
the same energy. However, the collapsed phase is sensi- 
tive to the short-distance form. That is, the system with 
repulsive short-distance force reenters the interval of pos- 
itive specific heat at a higher energy than those without 
such forces (see left panel of Fig. ||). 

Furthermore, the collapsed phase resulting from a sim- 
ulation with a Plummer potential is hotter and denser 
than those resulting from a simulation with repulsive small 
scale forces. This can be seen in Fig. [| as well. 

For a smaller softening length theory expects a phase 
transition at higher energies. Yet, already in simulations 
with a mild softening (e = 0.05) the collapsing transition 
occurs straight after the system has entered the interval 
of negative specific heat, which is shortly after the system 
has become self-gravitating. Thus, the collapsing transi- 
tion can not occur at significantly higher energies for a 
smaller softening length. 

However, for systems with smaller softening length, the 
energy at which the system reenters the zone of positive 
specific heat after the collapse changes, i.e., it is shifted 
to smaller energies. This is shown in Fig. ^, where two 
simulations with N = 160000 and e = 0.05 resp. e = 0.01 
are compared. 

Due to the higher particle number the dissipation time 
is in these simulations reduced to, Tdi S = 4rd yn . Thus, dy- 
namical equilibrium is not as perfect as previously. That 
is, the acceleration of the moment of inertia / attains a 
maximum value of, I/T = 0.14, where T is the kinetic 
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Fig. 4. Comparison of theoretical predictions and simulated systems. The simulations are carried out with N — 8000 
and the softening length is, e = 0.05. Consequently, the volume filling factor is, Vg = 1. Left: Inverse tempera- 
ture (3 versus energy e for models with softened potentials. The solid line indicates the theoretical result with the 
Follana potential, n — 10. The other curves depict the evolution of the simulated systems. Crosses: repulsive po- 
tential (£ = 2/3). Dots: Plummer potential (£ = 0). The circle indicates the initial state of the simulations. The 
range of negative specific heat corresponds to the range where the slope of /3(e) is positive. Middle: The dotted lines 
describe the evolution of the Lagrangian radii n for the simulation with the Plummer softening potential (£ = 0). 
Each curve depicts the radius of a sphere containing a certain mass fraction. The different mass fractions are: 
AM /M = {5%, 10%, 20%, . . . , 80%, 90%, 95%}. The solid line depicts the theoretical 95%-Lagrangian radius for the 
Follana potential. The dashed vertical line depict the moment when the system becomes fully self-gravitating. The 
dash-dotted and the dotted line mark the interval of negative specific heat. Right: Idem for the simulation with the 
short range repulsive force. 




Fig. 5. Same as Fig. [| for two simulations with N = 160000 and with a Plummer potential £ = 0.0. Dots: e = 0.05, 
V& = 20.0. Crosses: e = 0.01, V« = 0.16. 



energy. However, in general dynamical equilibrium is well 
approximated and we expect, due to the experience with 
simulations with varying particle number and dissipation 
time, that the effect of the temporal deviation from equi- 
librium does not affect qualitatively the results presented 
in Fig. ||. 



The volume filling factor of the systems with N = 
160000 and e = 0.05 is, Vg — 20, meaning that the mass 
resolution is a factor 2.7 larger than the force resolution. 
Thus the system is less granular than in the previous sim- 
ulations with N = 8000 and V s = 1 (see Fig. |). Yet, this 
does not affect the simulated quasi-equilibrium states and 
the deviation from theoretical predictions remains. 



5.1. Discussion 

Some predictions made by analytical models, such as the 
interval of negative specific heat, agree with findings re- 
sulting from iV-body models. Yet, there appear also dis- 
crepancies, such as the way the collapsing phase transition 
develops. 

At first sight the deviating results are surprising. Yet, 
the small scale physics is different in the two models which 
may account for the discrepancies. Indeed, thermostatis- 
tics assumes a smooth density distribution and is thus not 
able to account for two-body relaxation effects in granular 
media. However, a certain degree of granularity is a inher- 
ent property of TV-body systems. This may then lead to 
discrepancies, especially in the unstable interval of nega- 
tive specific heat, where phase transitions may be sensitive 
to small scale physics. 
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Whereas thermostatistics is too smooth to account 
for microscopic physics in granular self-gravitating media 
such as the interstellar gas, two-body relaxation is often 
too strong in iV-body system due to computational limita- 
tions. Especially in high force resolution simulations where 
the force resolution is larger than the mass resolution. We 
refer to this in Sect. 6.1 where we discuss, among other 
things, the effect of the granularity on long-range correla- 
tions appearing in the interval of negative specific heat. 

A further discrepancy between nature, analytical mod- 
els and iV-body models may be due to the entropy used 



in gravo-thermal statistics (Taruya & Sakagami 2001). 
Indeed, the entropy used in analytical models to find equi- 
librium states via the maximum entropy principle is the 
extensive Boltzmann-Gibbs entropy, that is in fact not ap- 
plicable for non-extensive self-gravitating systems. 

Generalized thermostatistics including non-extensivity 



are currently developed (Tsallis 198? ; Sumiyoshi 2001 
Latora et al 



2001; Leubner 2001 



These formalisms 

suggest that non-extensivity changes not all but some 



of the cla ssical thermodynamical results (Tsallis 1999; 
Boghosian 199£), which agrees with our findings. 

Because currently it is a priori not known which ther- 
mostatistical properties change in non-extensive systems 
and consistent theoretical tools are not available, analyti- 
cal results must be considered with caution. 

6. Results: II. Transient Long-Range Correlations 
in Gravitational Collapses 

6.1. Effect of Dissipation 

Properties of equilibrium states differing from theoretical 
predictions in the interval of negative specific heat were 
found (see above). Let us, now study the nonlinear struc- 
ture growth in this interval for an increased dissipation 
strength, i.e., outside of equilibrium, when Tdi s Si Tff . 

As previously, a relaxed, unperturbed, confined N- 
body sphere with e = 1 serves as initial state for the 
simulations. 

In order to dissipate the energy different dissipation 
schemes are applied. Before we discuss the appearance 
of long-range correlations in unstable dissipative systems 
let us briefly discuss the effect of the different dissipation 
schemes on the global system structure. 

The different applied dissipation schemes (see Sect. 
3.3) lead to collapsed phases with different global struc- 
tures. That is, they have different mass fractions con- 
tained in the core and the halo. A typical ordering is, 

^global > ^dyn > Aocal, where, Pglobal, 2? dyn and 2?l OC al 

are the density contrasts resulting from simulations with 
a global dissipation scheme, dynamical friction and a lo- 
cal dissipation scheme, respectively. The density contrast 
is, V = log(p cm / po), where p cm and po are the center of 
mass density and the peripheric density, respectively. This 
means that for a global dissipation scheme almost all the 
mass is concentrated in a dense core, whereas a local dis- 
sipation scheme can form a persistent "massive" halo. 



The evolution of the mass distribution resulting from 
simulations with the different dissipation schemes is shown 
in Fig. ||. The collapse of the inner shells takes in all sys- 
tems about the same time. Yet, the uncollapsed matter 
distributed in the halo is for the global dissipation scheme 
less than 2%, whereas it is 10% for the local dissipation 
scheme. 

Next, temporary long-range correlations that develop 
in unforced gravitating systems are presented in depen- 
dence of the different system parameters. 

A sufficiently strong global dissipation, i.e., T^i s Ss Tg, 
of a gravitating system leads during the nonlinear phase of 
the collapsing transition to fragmentation and long-range 
phase-space correlations, so that the index, <5(r), of the 
velocity-dispersion size relation, a oc r 5 ^ r > , becomes posi- 
tive. 

Fig. ^ shows the evolution of the velocity-dispersion- 
size relation, i.e., of S for different dissipation strengths. 
The relations result from 3 different simulations with 
global dissipation schemes and different dissipation 
strengths. The dissipation strength is given through the 
parameter a. Here a — 1.0, 5.0, 9.0. This corresponds to 
r dis = 2.0, 0.4,0.2r ff . 

Whereas the velocities remain uncorrelated in the sim- 
ulation with a = 1.0, S becomes temporary positive over 
the whole dynamical range in the simulations with the 
stronger dissipation. 

The velocity correlations start to develop at largest 
scales after the system has become self-gravitating. After 
the system has entered the interval of negative specific 
heat the correlation growth is accelerated. It attains a 
maximum and finally disappears when the collapse ends. 
The dynamical range over which S > during maximum 
correlation is ss 2 dex. 

Correlations at small scales straight above the soften- 
ing length are stronger for a stronger energy dissipation. 

The maximum correlation established in the negative 
specific heat interval persists for O(0. 1) ts and is charac- 
teristic for the applied dissipation strength and softening. 
For instance, the simulation with the strongest dissipation 
develops during 0.3 Tg a roughly constant S over a range of 
1.5 dex that resembles Larson's relation. Yet, correlations 
at small scales decay rapidly and the index d becomes even 
negative at intermediate scales. 

The end of the collapse and with it the disappearance 
of the correlated velocity structure is marked by a diverg- 
ing negative specific heat cy — * — oo. This is shown in Fig. 

s 

Systems with dynamical friction and local dissipation 
develop also velocity correlations. 

Fig. ^ shows the correlations resulting from simulations 
with a local dissipation scheme. Contrary to the global 
dissipation scheme, a strong local dissipation does not ex- 
tend the collapse of the whole system. Thus, local friction 
forces, that are strong enough to develop correlations lead 
also to a fast collapse and correlations persist accordingly 
a short time compared to simulations with global dissipa- 
tion. 
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Fig. 6. Evolution of the Lagrangian radii for a system with global dissipation scheme, dynamical friction and a local 
dissipation scheme, respectively. The Lagrangian radii depict, here and in the following figures, spheres containing the 
following mass fractions: AM/M = {5%, 10%, 20%, . . . , 80%, 90%, 98%}. 
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Fig. 7. Evolution of the index, 5(r), of the velocity-dispersion-size relation, a oc S(r), during the collapsing transition 
for three simulations with global dissipation scheme and different dissipation strength. The correlations result from 
simulations that were carried out with 160000 particles. The solid, the dashed and the dotted vertical lines indicate 
the scope of application of the simulation with a — 1.0, a — 5.0 and a = 9.0, respectively. The lower cutoff is given 
by the softening length, that is here e = 0.01, and the upper cutoff by 2Rgo. The time is indicated above each panel. 
The corresponding evolution of the Lagrangian radii and the specific heat are shown in Fig. @. 



The corresponding evolution of the Lagrangian radii 
and the interval of negative specific heat are shown in 
Fig. 0. 

In the dynamical friction scheme, energy dissipation 
depends, as in the global dissipation scheme, on the ab- 
solute particle velocity. Thus a strong dynamical friction 
extends the gravitational collapse and the "lifetime" of 
the correlations. Yet, the observed phase-space correla- 
tions are weaker, compared to those appearing in simula- 
tions with global and local dissipation, respectively, mean- 



ing that the reduced dissipation strength for fast particles, 
accordingly to the dynamical friction scheme, destroys the 
correlations. 

6.2. Effect of Short Distance Regularization 

So far the velocity correlations in dependence of the differ- 
ent dissipative factors were studied. Next the effect of dif- 
ferent short distance regularizations of the Newtonian po- 
tential, removing its singularity are checked. That is, sim- 
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Fig. 8. Evolution of the Lagrangian radii (top) and the specific heat (bottom) for three simulations with global 
dissipation and different dissipation strength. The dashed line marks the moment, when the system becomes fully 
self-gravitating. The dash-dotted line indicates the moment when the system enters the zone of negative specific heat. 
The dotted line marks the "end" of the collapse. The evolution of the corresponding phase-space correlations are shown 
in Fig. @. 
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Fig. 9. Evolution of the index of the velocity-dispersion-size relation for three simulations with a local dissipation 
scheme and different dissipation strengths. The particle number is, N = 160000, and the softening length is, e = 0.01. 
The solid, the dashed and the dotted vertical lines indicate the scope of application of the simulation with (3 — 1000, 
(3 = 2000 and (3 — 5000, respectively. The parameters determin ing th e dependence of the local dissipation on the radius 
and the velocity are, 77 = 1 and £ = 1, respectively (see Sect 3.3.l[ ). The evolution of the corresponding Lagrangian 
radii and the interval of negative specific heat are shown in Fig. hd. 



ulations with and without short distance repulsive forces 
and with different dissipation strengths are carried out. 
Energy is dissipated via the global dissipation scheme. 

In Fig. [ll] the velocity correlations resulting from three 
simulations with three different regularizations are com- 
pared. The regularizations are characterized by two pa- 
rameters. Namely, the softening length e and the pa- 
rameter £, that determines the strength of the repul- 
sive force. The softening length and £ of the three sim- 



ulations, compared in Fig. Ill], are, (e = 0.01, £ = 0.0), 
(e = 0.01, f = 2/3) and (e = 0.05, £ = 0.0), where £ = 0.0 
means that a Plummer potential is applied and £ > 1/3 
means that short distance repulsive forces are at work (see 
Sect. 3.5). The dissipation strength is for all three simu- 
lations the same, a = 10, and the collapsing time is con- 
sequently the same as well (see Fig. |l2"| ). 

During the first Tg long-range correlations resulting 
from the three simulations are identical. Then the S(r) 
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Fig. 10. Evolution of the Lagrangian radii and the specific heat for three simulations with local dissipation and 
different dissipation strength. The interval of negative specific heat is depicted by the dash-dotted and dotted line. 
The dashed line marks the moment when the whole system becomes self-gravitating. Fig. || shows the evolution of the 
corresponding long-range correlations. 



starts to separate. Indeed, the repulsive forces cause an in- 
crease of 6 at small scales, compared to the simulation with 
the same softening length, but without repulsive forces. 
Yet, large scales remain unaffected. 

The index S(r) resulting from the simulation with the 
large softening length is larger compared to the other two 
simulations after one Tg. Also, the 5(r) curve is flatter 
after that time. 

The large difference between the simulation with the 
large softening length and the corresponding simulation 
with the small softening length is astonishing, because 
there is a clear deviation even at large scales. Well, the 
two simulations were carried out with the same particle 
number, N = 160000. Consequently, the volume filling 
factors are different. Namely, Vs = 25 and Vs = 0.16, 
meaning that for the large softening length, i.e., the large 
Vg, the mass resolution is larger than the force resolution 
and vice versa for the small softening length. Whereas a 
small volume filling factor describes a granular phase, a 
large volume filling factor describes rather a fluid phase. 
Thus the deviation of the velocity correlations may mainly 
be due to the different volume filling factors and not due 
to the different softening lengths. 

In order to check this some complementary numerical 
experiments are carried out, whose results are presented 
subsequently. 

6.3. Granular and Fluid Phase 

The upper panels in Fig. |l3| show the velocity-dispersion- 
size relation, that develop in simulations with different vol- 
ume filling factors, but equal softening length, e = 0.031. 
The particle numbers are, TV = 6400, N = 32000 and 
N = 160000. Thus, the volume filling factors are, Vs = 
0.2, Vs = 1.0 and Vg = 5.0, respectively. Strongest phase- 
space correlations appear in the simulation with the high- 
est mass resolution, i.e., for the "fluid phase" , and weakest 
correlations appear in the "granular phase" . The flattest 
S(r) curve results from the simulation in which force and 
mass resolution are equal. 

This result suggests that the volume filling factor Vs 
(and not the softening length) is the crucial parameter 



determining the correlation strength in simulations with 
equal dissipation strength, as long as a substantial part of 
the system forces stem from interactions of particles with 
relative distances larger than the softening length. 

The lower panels in Fig. [l3] show, for the same simu- 
lations, the evolution of the index, D(r), of the mass-size 
relation, M oc r D ^'\ The deviation from homogeneity is 
strongest for the "granular phase" and weakest for the 
"fluid phase" . Thus the order of the correlation strength 
in space is inverse to the order of the correlation strength 
in phase-space. Such an inverse order is expected in self- 
gravitating system with |?7| oc T, where U is the potential 
energ y and T is the k inetic en ergy (Pfenniger & Combes 



1994| ; Pfenniger [199§ Combes |1999|) 



The corresponding Lagrangian radii and the interval 
of negative specific heat are shown in Fig. [l4|. 

6.4. Initial Noise 

Assuming a constant softening length, different volume 
filling factors, Vs, mean different particle numbers, N, 
that introduce different Poisson noise, ~ \/~N. Then, the 
statistical roughness of the initial uniform Poisson particle 
distribution decreases like ~ 1/y/N. 

Thus the initial roughness of the latter simulations (see 
Fig. [l^) differ from each other by a factor ~ 2.2 and one 
might suppose that the different correlation strengths in 
these simulations are the result of the different statistical, 
initial roughness. 

In order to check this possibility, long-range correla- 
tions are compared that result from two simulations with 
Vs — 1.0 and with statistical roughness that differs by a 
factor ~ 2.2. The results are presented in Fig. [l5| and the 
corresponding Lagrangian radii are shown in Fig. [l^. The 
simulation with N = 32000 was already presented in Fig. 
|l3| . It is now compared with a simulation with stronger 
initial roughness. 

Despite the different initial roughness, long-range cor- 
relations in phase-space are almost identical for the two 
simulations. 

As regards fragmentation, additionally to the parame- 
ter of the mass-size relation, the mass distribution in space 
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Fig. 11. Evolution of the velocity correlations for three simulations with global dissipation scheme and different 
potential regularizations. That is, the softening length and the form of the softening, respectively, change from 
simulation to simulation. For £ = 0.0 the particle potential is a Plummer. £ = 2/3 means that there are short distance 
repulsive forces at work. The particle number is, N=160000. The solid, the dashed and the dotted vertical lines 
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respectively. The evolution of the corresponding Lagrangian radii are shown in Fig. H3. 
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Fig. 12. Same as Fig. [To] for three simulations with different potential regularizations, i.e., different softening 
lengths and forms of the softening potential, respectively. The evolution of the corresponding phase-space correlations 
are shown in Fig. O. 



were compared. We find actually a stronger fragmentation 
for small particle numbers. Yet, the Vg seems to be the 
crucial parameter for the fragmentation strength in Fig. 

Consequently, the different correlation strengths in 
Fig. [l3| can not be accounted for by Poisson noise, but are 



mainly the result of the different volume filling factors, or 
more precisely, due to the different ratio of force and mass 
resolution. Thus dark matter clustering in high-resolution 
cosmological TV-body simulations in which the force res- 
olution is typically an order of magnitude smaller than 
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Fig. 13. Upper panels: Evolution of the index, of the velocity-dispersion-size relation, a oc r 5 ( r ), resulting from 
three simulations with different volume filling factor that describe a granular phase, Vg = 0.2, a fluid phase, 
Vs = 5.0, and a intermediate state, Vg = 1.0. The particle numbers are iV = 6400, N = 160000 and N = 32000, 
respectively. The softening length is e = 0.031. The solid, the dashed and the dotted vertical lines indicate the scope of 
application of the simulation with Vg — 0.2, Vg — 1.0 and Vg = 5.0, respectively. The scope is give by, e < r < 2Rqq. 
The corresponding Lagrangian radii are presented in Fig. [l4|. Lower panels: The evolution of the index, D(r), of the 
mass-size relation, M oc r Dl ^ r \ The solid, the dashed and the dotted vertical lines indicate the scope of application of 
the simulation with Vs — 0.2, Vs — 1.0 and Vg = 5.0, respectively. The scope is give by, e < r < Rgo/2. 



the mass resolution may be too strong compared with the 
physics of the system (Hamana et al. 2001). 



6.5. Free Fall of Cold Dissipationless Systems 

Here, long-range correlations are discussed that develop in 
cold, gravitationally unstable systems without energy dis- 
sipation. The correlations strength appearing in such dissi- 
pationless systems depends on the initial ratio between ki- 
netic and potential energy, a = T/\U\, i.e., on the number 
of thermal Jeans masses given through M/Mj — 2a~ 3 / 2 . 

This is shown in Fig. |l7|. For a = 0.01 the index 8 of 
the velocity-dispersion-size relation remains zero at small 
scales during the whole free fall, whereas for a = 0.0 the 



index becomes, 8 > 0, over the whole dynamical range, 
during a period of ~ 0.7rgf. That is, 350 Mj are not suffi- 
cient to develop small scale phase-space correlations in a 
dissipationless system. 

In order to show the dependence of the result on the 
initial Poisson noise, the absolutely cold simulation with 
N = 160000 is compared with a N = 32000-body simula- 
tion. As above, it follows that the initial roughness of the 
two systems differs by a factor ~ 2.2. 

Fig. [l7| shows that the non-equilibrium structures re- 
sulting form the simulations with equal a but unequal par- 
ticle number differ from each other during the whole free 
fall. Indeed, during the first half of the free fall time the 
two initial conditions produce velocity correlations that 
differ on small and large scales. After t = 0.5 Tg the dif- 



14 



D. Huber and D. Pfenniger: Correlations in Self-Gravitating Systems 




ferences approximately disappear. However, the behavior 
of the spatial correlations (see lower panels of Fig. [l7]) is 
inverse, meaning that differing spatial correlations appear 
after 0.5 Tg and persist for the rest of the free fall. 

These results suggest that non-equilibrium structures 
appearing in cold systems or in systems with very effective 
energy dissipation depend more strongly on initial noise 
than those appearing in warm systems with less effective 
dissipation (see also Fig. |l5| ). 

The evolution of the Lagrangian radii during the free 
fall of the cold, dissipationless systems are shown in 
Fig. " 
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Because these systems are adiabatic and self- 
gravitating during the whole simulation, the correspond- 
ing intervals are not plotted in the figure. 



6.6. Discussion 

Long-range spatial and phase-space correlations appear 
naturally during the collapsing phase transition in the in- 
terval of negative specific heat if the energy dissipation 
time is, Tdi S ts, so that the time-scale of correlation 
growth is smaller than the time-scale of chaotic mixing 
which is, ~ 1/A o c Tg, where A is the maximum Liapunov 
exponent (Miller |l994j) . 

Actually, the details of the long-range correlations de- 
pend on the applied dissipation scheme, but there are also 
some generic properties. That is, phase space correlations 
start to grow at large scales, whereas spatial correlations 
seem to grow bottom-up. Moreover, there is an upper limit 
for the index of the velocity-dispersion-size relation within 
the dynamical range, namely, 1, 

Besides, the dissipation strength, and initial condi- 
tions, the volume filling factor Vg is a crucial parameter 
for the correlation strength. That is, phase-space correla- 
tions are strong in the fluid limit and weak for a granular 
phase. The behavior of the spatial correlations is exactly 
the other way around. 

The softening length does not affect correlations within 
the dynamical range. Yet, sub-resolution repulsive forces 
affect correlations on small scales above the resolution 
scale. Of course, this does not hold for the onset of the 
correlation growth, but only after sufficient particles have 
attained sub- resolution distances. 



The above considerations suggest that the found cor- 
relations are physically relevant and not a numerical arti- 
fact. 

The long-range spatial and phase-space correlations 
appearing during the collapsing transition are qualita- 
tively similar to the mass-size and the velocity-dispersion- 
size relation observed in the ISM (e.g. Blitz & Williams 



1999; Chappell & Scalo 2001; Fuller & Myers 1992), show- 



ing that in the models gravity alone can account for ISM- 
like correlations. 

Furthermore, the time-scale of the correlation lifetime 
is the free fall time, Tg, which is consistent with a dy- 
namical scenario in which ISM-structures are highly tran- 
sient (Vazquez-Semadeni 2002; Larson 2001; Klessen et 
al. 200C), which is related with rapid star formation and 



short molecular cloud lifetimes (Elmegreen 2000), that is, 
the corresponding time-scales are about an order of mag- 
nitude smaller than in the classical Blitz & Shu (1980) 
picture. 

7. Results: III. Permanent Energy-Flow 

7.1. Nonequilibrium Structures in Systems Subject to 
an Energy Flow 

The numerical experiments presented above showed that 
dissipative self-gravitating systems fragment and establish 
long-range correlations outside of equilibrium in the inter- 
val of negative specific heat. These transient correlations 
persist for 1 — 2 Tg times. 

Here we check if self-gravitating systems can estab- 
lish persistent long-range correlations if they are main- 
tained continuously outside of equilibrium by a permanent 
energy-flow. That is, the dissipated energy is continuously 
replenished by time-dependent potential perturbations. 

Both simulations of granular and fluid phases are car- 
ried out. A typical parameter set, describing a granular 
system is, N = 10000, e = 0.0046, Vg = 0.001 and 
N = 160000, e = 0.0315, Vg = 5 are here typical for a 
fluid phase. Yet, the parameters are not fixed and the ef- 
fect on the evolution is studied when parameters change. 
Parameters, controlling energy-flow and interaction poten- 
tial, as well as their ranges are indicated in Table |l[ 

The applied potential perturbations imitate massive 
objects passing in the vicinity on time-scales, T port Td yn - 
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Fig. 15. Upper panels: Evolution of the phase-space correlations resulting, resulting from two simulations with identical 
volume filling factor, Vff = 1.0, but unequal statistical, initial roughness, ~ 1/^/N. The particle number and the 
softening length of the first simulation (dots) arc N = 32000 and e = 0.031, respectively. The second simulation 
(triangles) is carried out with N = 6400 and e = 0.054. The solid and the dashed lines indicate the scope of application 
given by e < r < 2Rgo- The evolution of the corresponding Lagrangian radii is shown in Fig. [l6| Lower panels: The 
evolution of the spatial correlations. The scopes of application, depicted by the vertical lines, is, e < r < Rqq/2. 



The perturbations induce primarily ordered particle mo- 
tions. Then gravitational interactions lead to a conversion 
of the bulk kinetic energy to random thermal motion. The 
energy injection due to such a forcing scheme can be quite 
regular until a plateau is reached. 

Energy injection prevents a system from collapsing and 
maintains an approximately statistically steady state dur- 
ing ^5 — 15 Tdyn when energy dissipation is balanced ap- 
propriately by large scale potential perturbations. These 
states do not feature any persistent long-range correla- 
tions. Yet, they develop a temperature structure that is 
characteristic for the applied dissipation scheme. This is 
shown in Fig. [l^, where the evolution of two granular sys- 
tems subjected to an energy-flow is presented. One system 
dissipates its energy by a global dissipation scheme (top), 
the other by a local dissipation scheme (bottom). 



The system with the global dissipation scheme is 
nearly thermalized during almost the whole simulation, 
whereas the local dissipation scheme leads to a permanent 
positive temperature gradient, that is inverse compared to 
stars and resembles those of the ISM where dense, cool 
mass condensation are embedded in hotter shells. 

If dissipation dominates energy injection the system 
undergoes in general a mono-collapse, that is, a collapsed 
structure is formed in which a part of the system mass 
is concentrated in a single dense core and the rest is dis- 
tributed in a diffuse halo. However, systems with a rather 
fluid phase may develop several dense cores moving in a 
diffuse halo, when they are subjected to an appropriate 
energy flow. In the course of time the number of clumps 
varies, but a non-mono-clump structure persists for some 
Tff (see Fig. ||c]). These systems may even develop persis- 
tent phase-space correlations (see Fig. ^l| ). 
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Fig. 16. Same as Fig. [LO] for two simulations with identical volume filling factor and unequal initial roughness, 
~ 1/Vn. The evolution of the corresponding long-range correlations is shown in Fig. |l5| 



However, the clumps result not from a hierarchical 
fragmentation process, but they are formed sequentially 
on the free fall time scale. Furthermore, the clumps are 
so dense that their evolution is strongly influenced by the 
applied regularization. That is, the evolution over several 
Tff depend on the numerical model that does not represent 
accurately small-scale physics. 

In order to impede gravitational runaway that may 
hinder the formation of complex non-equilibrium struc- 
tures within the given dynamical range, different measures 
are taken, such as short distance repulsion and the appli- 
cation of the dynamical friction scheme, where the friction 
force F — > for v — > oo. 

Yet, despite these measure, we only find either nearly 
homogeneous structures in systems where energy injection 
prevails, or systems dominated by unphysical clumps in 
case of prevailing energy dissipation. 

Up to now, mainly the effect of several different dis- 
sipation schemes were discussed. However, the effect of a 
modified forcing scheme is also checked. That is, a power- 
law forcing scheme is applied, which injects energy at dif- 
ferent frequencies. This forcing scheme is a modification 
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of those presented in Sect. 3.4 and reads, $ pGrt (w) oc w", 
where v = — 4, . . . , 4. Yet, also such a forcing scheme can 
not induce a phase transition to complex non-equilibrium 
structures and the resulting mass distribution corresponds 
to the those described above. 



7.2. Discussion 

Actually the "clumpy" structure in Fig. and the cor- 
responding phase-space correlations shown in Fig. ^l] do 
not represent real physics, nevertheless, they show that it 
is principally possible to maintain spatial non-equilibrium 
structures and long-range phase-space correlations in a 
perturbed, dissipative, self-gravitating system over sev- 
eral dynamical times. Thus, it can not be excluded that 
in the future, with a better representation of microscopic 
physics and forcing mechanisms at work in the ISM, mod- 
els including self-gravity may produce complex nonhomo- 
geneous structures in a statistical equilibrium, i.e., persis- 
tent patterns formed by transient structures. 

However, at present models of dissipative, self- 
gravitating systems, cannot produce such structures on 



On larger scales, gravity gives rise to persistent non- 
equilibrium structures in cosmological and shearing box 
simulations. A common denominator of these models is 
that their time dependent boundary conditions are given 
by a scale-free spatial flow counteracting gravity. Let us 
discuss this more precisely. 

In cosmological and shearing box models time- 
dependent boundary conditions create relative particle ve- 
locities that are inverse to gravitational acceleration and 
increase with particle distance, v oc r. In the shearing 
box model, the relative azimuthal particle velocity due to 
the shear flow is, vg oc r c , where r c is the radial parti- 
cle distance in cylinder coordinates. In cosmological mod- 
els, the relative particle velocity induced by the Hubble 
flow is, v r oc r, where r is the relative particle distance in 
Cartesian coordinates. 

These relations and consequently the corresponding 
flows are scale-free. The fact that the shear flow affects 
only the azimuthal velocity component may then account 
for the characteristic spiral arm like structures found in 
shearing box experiments, differing from those found in 
cosmological models, where the isotropic Hubble flow gives 
rise to in average isotropic non-equilibrium structures. 

The models studied in this paper are not subject to 
a scale-free spatial flow counteracting gravity and persis- 
tent long-range correlations of astrophysical relevance do 
not appear. This may suggest that in situations where 
gravitational runaway is allowed, matter that has passed 
through a collapsing transition, has to be replenished at 
large scales in order to attain a statistical equilibrium state 
of transient fragmentation. 

8. Conclusions 

First, equilibrium states of iV-body models were compared 
with analytical models. Subsequently, the findings result- 
ing from this comparison are summarized: 

— On the one hand, equilibrium properties of ./V-body 
models agree with predictions made by analytical mod- 
els. An example is the energy interval of negative 
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Fig. 17. Upper panels: Phase-space correlations that develop during the free fall of three cold, dissipationless systems. 
All systems have the same volume filling factor, Vg — 0.016. The parameter a indicates the ratio between kinetic and 
potential energy. The solid, the dashed and the dotted lines indicate the scope of application, e < r < 2i?go, of the 
simulation with (a = 0.01, N = 160000), (a = 0.0, N = 160000) and (a = 0.0, N = 32000), respectively. The evolution 
of the corresponding Lagrangian radii is shown in Fig. Lower panels: Evolution of the spatial correlations. The 
scopes of application, depicted by the vertical lines, is, e < r < Rqq/2. 



specific heat. One the other hand, discrepancies were 
found, such as the way the collapsing phase tran- 
sition, separating a high-energy homogeneous phase 
from a low-energy collapsed phase, develops in the in- 
terval of negative specific heat. These discrepancies 
suggest: 1.) Small scale physics becomes relevant for 
the system evolution when the growth of singularities 
triggered by gravitational instabilities is allowed. 2.) 
Analytical models based on the Gibbs-Boltzmann en- 
tropy are not strictly applicable to non-extensive self- 
gravitating systems. Yet, not all of the equilibrium 
properties found by maximizing the Gibbs-Boltzmann 
entropy are expected to change if a fully consistent, 
generalized thermostatistical theory is applied. 

Second, the collapsing transition was studied in systems 
with strong dissipation. The findings are: 



— Dissipative self-gravitating systems develop outside of 
equilibrium, in the interval of negative specific heat, 
transient long-range correlations. That is, fragmen- 
tation and nonequilibrium velocity-dispersion-size re- 
lations, with striking resemblance to those observed 
in the ISM, appear during the collapsing transition, 
when the dissipation time is shorter than the dynami- 
cal time. This suggests that nonequilibrium structures 
in self-gravitating interstellar gas are dynamical and 
highly transient. 

— Besides the dissipation strength and the initial noise, 
the granularity turns out to be a crucial parameter for 
the strength of the resulting long-range correlations, 
substantiating the importance of a coherent mass and 
force resolution. That is, phase-space correlations are 
stronger in the fluid limit than in a granular phase. The 
opposite holds for spatial correlations. This substanti- 
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Fig. 18. Same as Fig. [K] for the free fall of three cold, dissipationless systems with identical volume filling factor 
but unequal a and initial roughness, respectively. The evolution of the corresponding long-range correlations is shown 
in Fig. §3. 




t„ t„ 5 10 15 5 10 15 



Fig. 19. The two left panels: Evolution of the Lagrangian radii (left) and the Lagrangian dimensionless temperature 
(middle left) of a system subjected to an energy-flow and with global dissipation scheme. The curves depict certain 
mass-fractions, and its temperatures, respectively. The mass fractions, contained in spheres, centered at the center of 
mass, are: AM/M = {5%, 10%, 20%, . . . , 80%, 90%, 98%}. The triangles (left panel) indicate the temperature of a mass 
fraction of 5% and the crosses corresponds to a mass fraction of 98%. The simulation is carried out with N = 10000 
and e = 0.0046. Consequently, Vg = 0.001, meaning that the system is highly granular. The dynamical range is, 2.3 
dex. The two right panels: Same as above for a system subjected to an energy-flow and with a local dissipation scheme. 



ates The inverse behavior of fragmentation strength 
and phase-space correlation strength is found in all 
simulations and is typical for self-gravitating systems. 

Finally, systems subject to a permanent energy-flow were 
studied. We find: 

— Typically driven dissipative systems evolve to a 
high-energy homogeneous phase or undergo a mono- 
collapse. Yet, model systems with a local energy dis- 
sipation can develop persistent phase-space correla- 
tions, but a persistent, hierarchical fragmented struc- 
ture is not observed. This suggests, that matter that 
has passed through a collapsing transition has to be 
replenished at large scales in order to maintain a hier- 
archical structure at molecular cloud scales. 
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